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We present a comprehensive study of TeV black hole events in Earth's atmosphere originated by 
cosmic rays of very high energy. An advanced fortran Monte Carlo code is developed and used to 
simulate black hole extensive air showers from ultrahigh-energy neutrino-nucleon interactions. We 
\ investigate the characteristics of these events, compare the black hole air showers to standard model 

f""^ . air showers, and test different theoretical and phenomenological models of black hole formation and 

' evolution. The main features of black hole air showers are found to be independent of the model 

I considered. No significant differences between models are likely to be observed at fluorescence 

telescopes and/or ground arrays. We also discuss the tau "double bang" signature in black hole air 
showers. We find that the energy deposited in the second bang is too small to produce a detectable 
peak. Our results show that the theory of TeV-scale black holes in ultrahigh-energy cosmic rays 
leads to robust predictions, but the fine prints of new physics are hardly to be investigated through 
atmospheric black hole events in the near future. 
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I. INTRODUCTION 



The study of super-Planckian coUisions dates back to the late 80's (ij. Today's renewed interest stems from 
the possibility that the fundamental scale of gravity may be much lower than the observed gravitational scale Q . In 
i-S^ ' braneworld scenarios, the observed weakness of the gravitational field is due to the "leakage" of gravity in the extra 
dimensions: Standard model (SM) fields are constrained in a four-dimensional submanifold, whereas gravitons are 
O i' allowed to freely propagate in the higher-dimensional spacetime Q. If the gravitational coupling constant is of the 
D , order of few TeVs, the physics of super-Planckian collisions could soon be detected through observation of subnuclear 
black holes (BHs) and other extended objects, such as branes, in particle colliders ^ ^ or ultrahigh-energy cosmic 
j> \ ray (UHECR) observatories 0jIEHj0|- (For reviews and extra references, see Refs. Pllllll). 

. The semiclassical limit of super-Planckian scattering suggests that the cross section for creation of a BH or brane 
' with radius R is approximately given by the geometrical black disk UBD{s,n) = TrR'^{s,n), where ^/s is the center 
5^ , of mass (cm.) energy of the colliding quanta and n is the number of extra dimensions. Gravitational objects with 
mass of order of the fundamental gravitational scale M^, have radius of order M~^. In symmetric compactification 
models, the size of extra dimensions is much larger than . (For conventions, see Ref. Thus the spherical 

approximation is justified; the geometry of nonperturbative objects is that of a n-dimensional BH. The spherical 
approximation breaks down for asymmetric compactifications, where some of the extra dimensions have size of order 
of the fundamental Planck scale. In that case, the geometry of nonperturbative objects is that of strings and branes 

UHECRs arc attractive because of their high cm. energy. The nucleon-nucleon cross section for formation of BHs 
and branes is very small compared to other SM hadronic processes. The neutrino-nucleon cross section for BH or brane 
formation may be higher than the cross section of the SM process, thereby giving interest to neutrino interaction. 
Under the most favorable circumstances, the cross section for BH formation at the TeV scale reaches millions of 
pb for neutrino-nucleon collisions in the atmosphere. The cross section for brane production is expected to be even 
larger. These results have led to the claim that UHECR detectors might observe BHs and probe Planckian physics. 
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Event observables would be the secondary products of the BH or brane after their formation, i.e. extensive air showers 
originated by field emission in the decay phase. Hawking evaporation provides an emission mechanism for BHs |l4j . 
SM quanta are emitted in the visible three-brane and can be detected. Branes may or may not evaporate, depending 
on their properties. However, the decay spectrum of massive excitations in string theories has been shown to be 
thermal [l5|. This suggests that BH and brane decay signatures may be similar. 

The observational signatures of BH events in the semiclassical approximation have been investigated in a number 
of recent publications. In Ref. the authors found that BH interactions generate different air showers from SM 
interactions. BH air showers tend to rise faster and have larger muon content. A BH event produces a hadronic air 
shower occurring at a much greater depth in the atmosphere, i.e., a very deeply penetrating hadronic air shower. 
However, the inability of realistic detectors to observe the first interaction point hides most of the difference between 
BH and SM air showers. Given that present observatories are not large enough to study a large number of neutrino 
events, discrimination of BH and SM events is likely not to be achieved in the near future. Another characteristic of 
BH air showers is r generation. Although the rate for these events is low if the ultrahigh-energy neutrino flux is at 
the level of the expected cosmogenic neutrinos, unusual "double-bang" air showers could signal a departure from SM 
interactions. 

Whereas the semiclassical picture seems reasonable, the actual physics of subnuclear BH formation could be very 
different. In the last year or two, significant advances in the understanding of microscopic BH formation and air 
shower evolution have appeared in the literature. It is thus timely and worthwile to re-examine the observational 
signatures of BH air showers. To this purpose, we developed a thorough fortran Monte Carlo (MC) code to simulate 
the air showers induced by BH formation in neutrino-air collisions, which includes these theoretical refinements 
The MC has the same structure of the MC used in Ref. 9]. The code generates observable secondaries from BH 
evaporation using the PYTHIA generator [l^. These secondaries are then injected into the AIRES simulator jl^ as 
primaries for the final air shower. 

The purpose of our study is threefold. Firstly, we want to confirm the main findings of the previous investigation. 
Secondly, we want to test various proposals of BH vs. SM air shower discrimination that have appeared in the 
literature, such as the r "double bang" effect Thirdly, we want to look for new ways of discriminating between 
different models of BH formation and evolution. 

Our analysis will show that the main characteristics of BH air showers are essentially independent of the details of 
BH evolution. Because of large uncertainties and statistical fluctuations in air shower detection, it is also practically 
impossible to discriminate between alternative models of BH formation and evaporation. For instance, we will show 
that there is no significant observational difference between a model of BH formation based on the semiclassical black 
disk and the trapped-surface model [l^ , or between a model of BH evaporation with final explosive decay and stable 
remnant These results limit significantly the use of BH air showers (if they exist) to probe details of "new 

physics" . We will also show that newly proposed signatures do not help in the task of discriminating BH vs. SM air 
shower detection. No observational trace of the "double bang" signature can be extracted from a realistic detector in 
the near future. 



II. BASICS OF BH FORMATION AND EVOLUTION 

In this section we briefly review the basics of BH formation and evolution, focusing on recent theoretical advances 
that have been included in the MC code. 



A. BH formation and cross section at parton level 

Thome's hoop conjecture states that a horizon forms when a mass M is compacted into a region with circum- 
ference smaller than twice the Schwarzschild radius R{M) in any direction. At subnuclear level, this can be achieved 
by scattering two partons (ij) on the brane with cm. energy ^/Sij > M and impact parameter b < R{M). This event 
can formally be described by the process ij — + BH -I- E{X), where E{X) denotes collisional energy that does not 
contribute to the BH mass. This energy includes a bulk component of gravitational radiation and perhaps non-SM 
gauge fields, and a brane component of SM fields. If E(X) is zero, the hoop conjecture implies that the cross section 
for BH production is independent of the impact parameter (as long as 6 < R{M)) and equal to the geometrical black 
disk aBoisij, n). If E{X) ^ 0, the cross section depends on the impact parameter, and is expected to be smaller than 
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the black disk cross section. It is worth stressing that this picture is correct only if the BH is larger than the Compton 
length of the colliding quanta. (For discussions on the effect of wave packet size on the BH formation process, see 
Ref. [iJl-) A precise calculation of the coUisional energy loss is essential to understand BH formation. 

Many papers have been devoted to improve or disprove the hoop conjecture. The most popular model is currently 
the trapped-surface model although alternative techniques have been explored 13]. The trapped-surface 

approach gives an upper bound on the gravitational component of E{X) by modelling the incoming partons as two 
Aichelburg-Sexl shock waves The Aichelburg-Sexl wave is obtained by boosting the Schwarzschild solution to the 
speed of light at fixed energy. The resulting metric describes a plane-fronted gravitational shock wave corresponding to 
the Lorentz-contracted longitudinal gravitational field. The parton scattering is simulated by superposing two shock 
waves travelling in opposite directions. The union of these shock waves defines a closed trapped-surface that allows 
to set a lower bound on the BH mass. The collisional energy loss depends on the impact parameter and increases 
as the number of spacetime dimensions increases. The BH mass monotonically decreases with the impact parameter 
from a maximum of about 60-70% of the cm. energy for head-on collisions. 

The trapped-surface result is consistent within one order of magnitude with the hoop conjecture. However, the 
partons are assumed to be pointlike, massless, spinless, and electrically neutral. The pointlike assumption fails for 
directions transversal to the motion p^ . Colliding partons generally have spin and charge. While size and spin 
effects are expected to be mostly relevant around the Planck energy, charge effects could dominate at higher energy. 
It should also be kept in mind that the trapped-surface model provides only a lower bound on the BH mass. An 
accurate estimate of the gravitational collisional energy loss would require the use of the full non-linear Einstein 
equations in higher dimensions. Since this is a virtually impossible task, alternative approximated models have been 
investigated. The gravitational energy emission in a hard instantaneous collision can be computed in the linearized 
limit 27]. This approach suggests that the trapped surface method overestimates the gravitational energy emitted in 
the process. For head-on collisions, the instantaneous method predicts the gravitational energy loss to be only about 
10% of the cm. energy. This result is in agreement with a perturbative calculation modelling the parton-parton 
collision as a plunge of a relativistic test particle into a BH with mass equal to the cm. energy |28| . 

In conclusion, known methods for the estimate of the gravitational loss in relativistic scattering at parton level give 
a BH mass ranging between 60% and 100% of the cm. energy. Today, the trapped-surface value and the black disk 
value can be considered as the lower and upper bounds on the BH mass, respectively. 



B. Cross section at nucleon level 



The total cross section for a super-Planckian event involving a nucleon is obtained by integrating the above cross 
section over the parton distribution functions. BHs formed in a neutrino- nucleon collision may dominate over the 
SM processes and stand a fair chance of detection. On the contrary, the branching ratio of the BH cross section in a 
nucleon- nucleon collision is ~ 10~^. Therefore, BH detection in nucleon- nucleon interactions cannot be achieved with 
current and next generation detectors due to the low flux of UHECRs. 

If the BH mass depends on the impact parameter, the generally accepted formula for the total cross section of the 
neutrino-nuclcon process is 

UyN-yBH I 2zdz / dxqi{x,~Q'^) FaBD{xs,n) , (1) 

where qi{x, —Q^) are the Parton Distribution Functions (PDFs) [2^l30l ] with four-momentum transfer squared —Q^, 
and fraction of the nucleon's momentum carried by the i-th parton ^/x. z is the impact parameter normalized to its 
maximum value and Xm — M'^^^j {sy'^{z)), where y{z) and Mmin are the fraction of cm. energy trapped into the BH 
and the minimum-allowed mass of the gravitational object, respectively. F is a form factor. The total cross section 
for the black disk model is obtained by setting F = 1 and y'^{z) = 1. 

Different sets of PDFs are defined in the literature. The PDFs are not known at energies above the TeV and for 
values of momentum transfer expected in BH formation. Equation is usually calculated by imposing a cut-off at 
these values. The PDFs also suffer from uncertainties at any momentum transfer (~ 10%) |^ and from the ambiguity 
in the definition of Q jsj • The momentum transfer is usually set to the BH mass or the inverse of the Schwarzschild 
radius. Although recent literature inclines toward the latter, there are no definite arguments to prefer either one or 
to exclude alternative choices. The uncertainty due to the ambiguity in the definition of the momentum transfer is 
about - 10-20% 0. 
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The form factor and the amount of trapped energy depend in principle on energy, gravitational scale, geometry and 
physical properties of the extra dimensions and gravitational object. The trapped-surface method gives numerical 
values of order unity for these quantities. (See Refs. 1^ and discussion above). However, these results depend 
on the way the trapped-surface is identified. Other models 32^ give values which are more or less consistent with the 
trapped-surface method. With the lack of further insight, it is common practise in the literature to either choose the 
trapped-surface result or the simple black disk model. 

The lower cutoff on the fraction of the nucleon momentum carried by the partons is set by the minimum-allowed 
(formation) mass of the gravitational object, Mmm- This threshold is expected to roughly coincide with the mass 
for which the semiclassical description is valid. This conclusion is motivated for spherically symmetric BHs by the 
following argument |33|: For Mmin/M^ > few, the Hawking entropy of the BH should be large enough to neglect 
strong gravitational effects. The semiclassical results are then extrapolated for smaller values with the assumption that 
the BH or its Planckian progenitor decays on the brane. However, this argument is based on Hawking's semiclassical 
theory and may not be valid at energies equal to few times the Planck mass. For example, the existence of a minimum 
spacetime length l„i implies the lower bound on the BH mass JJt (^3 • 

Mmi = (2V^/„./Af.)"+' M. . (2) 

BHs with mass less than M^i do not exist, since their horizon radius would fall below the minimum-allowed length. 
At fixed M^,, the minimum- allowed mass grows as a power of l^^. For n = 6 or 7 and Im^I* ^ 1, it follows 



C. BH evolution 



It is believed that the decay of microscopic BHs happens in four distinct stages: I. radiation of excess multipole 
moments (balding phase); II. spin-down; HI. Hawking evaporation; IV. final explosion or formation of a BH remnant. 
Phases Till rely on semiclassical results, provided that the entropy is sufficiently large. Phase IV is in the realm of 
quantum gravity. 

Although some progress has been made, the understanding of balding phase and spin-down phase is still fragmentary. 
For example, the emission of radiation from a (n -I- 4)-dimensional rotating BH on the brane is known only for spin-0 
fields ^5|. Due to these limitations, some of these theoretical results cannot be implemented in MC simulations at the 
present stage. Moreover, phase I is not expected to lead to a significant amount of energy loss and the observational 
uncertainties (see below) are likely to dominate the theoretical uncertainties in phase II. 

Many papers have been devoted to the investigation of the Hawking phase. Although several analytical and 
numerical results have been obtained |36j| . from the viewpoint of numerical simulations the situation is similar to the 
balding and spin-down phases. For instance, greybody factors for the graviton are not fully known in (n-|-4)-dimensions 
even for the spherically symmetric BH. This precludes their use in numerical codes, where a consistent use of greybody 
factors is required. The field content at trans-Planckian energies is also not known. Onset of supersymmetry, for 
example, could lead to other evaporation channels for the BH and large emission of non-SM or undetectable quanta 
during the decay phase. Finally, quantum effects may also affect the emission of visible quanta on the brane. 

Quantum corrections to the Hawking phase can be phenomenologically described by assuming the existence of a 
minimum length of the order of the Planck length [37j| . The existence of a minimum scale is a common consequence 
of most (if not all) theories of quantum gravity such as string theory, non-commutative geometry and loop quantum 
gravity. The presence of a cutoff at small spacetime distances leads to a modification of the uncertainty principle at 
Planck scales. Since the Hawking thermodynamical quantities can be derived by applying the uncertaint y p rinciple 
to the BH, the existence of a minimum length leads to corrections in the thermodynamical quantities |l9ll34| . 

At the end of the Hawking phase, the BH is expected to either non-thermally decay in a number Up of hard quanta 
or leave a remnant. In either case we must content ourselves with a phenomenological description, due to the lack of 
a theory of quantum gravity. The final non-thermal decay is usually described by setting a cutoff on the BH mass of 
the order of the Planck mass, Qmin ~ M-i,, and then equally distributing the energy Qmm to a number Hp of quanta. 
Since the decay is non-thermal, and in absence of any guidance from a theory of quantum gravity, the quanta are 
democratically chosen among the SM degrees of freedom. Note that Qmin does not necessarily coincide with Mmin- 
The former gives the threshold for the onset of quantum gravity effects, whereas the latter gives the minimum- allowed 
mass of the classical object. From the above definitions, it follows Mmin ^ Qmin- The existence of a minimum length 
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gives a natural means to set Qmin- In that case, the modified thermodynamical quantities determine the endpoint of 
Hawking evaporation when the mass of the BH reaches Mj^i. This mass can be identified with the mass of the BH 
remnant [Tll33|. 



III. BH GENERATOR 

In this section we fist the main characteristics of the MC generator used in the simulations |0| . The physics of BH 
formation and decay is determined by the following set of external parameters and switches in the MC code: 

1. Fundamental Planck scale; 

2. Number of extra dimensions; 

3. Gravitational loss at BH formation and gravitational loss model; 

4. Minimum BH mass at formation; 

5. Quantum BH mass threshold at evaporation; 

6. Number of final quanta at the end of BH decay; 

7. Momentum transfer model in parton collision; 

8. Conservation of electromagnetic (EM) charge; 

9. Minimum spacetime length. 

The above parameters are briefly explained below. A more detailed explanation can be found at the MC generator 
web site Hal. 



A. BH formation and parton cross section 

The MC does not require any lower or upper bound on the Planck mass M*. However, experimental constraints 
exclude values of A/^ < 1 TeV and large values of M^, do not allow BH formation in the atmosphere. Therefore, Mi, 
must be chosen with caution. Since n = I and n = 2 are excluded experimentally, and most of the theoretical models 
are limited to n < 7, the number of extra dimensions n ranges from 3 to 7. 

The MC includes three models for BH formation and cross section: Black disk, Yoshino-Nambu (YN) trapped- 
surface model jl^, and Yoshino-Rychkov (YR) improved trapped-surface model |^. This allows a comparison 
between air showers based on the black disk model [Sj and air showers generated by BHs with significant gravitational 
loss at formation. Observable differences between different models of BH formation can be investigated, as suggested 
in Ref. |2|. 

The minimum BH mass Mmin is set in units of Mi, or Mmi (if a minimum length is present, see below). This 
parameter is always larger than one, i.e. Mmin > Max(M+, Mmi)- 



B. Total and differential cross section 

The distribution of the initial BH masses is given by the differential cross section da / cLMbhi where Mbh — \/xs- 
The MC uses the (stable) cteqS PDF distribution [30l l38j | . Since the use of different PDF distributions produces an 
insignificant uncertainty in the total and differential cross sections, other PDF distributions are not implemented in 
the MC. The uncertainty due to the choice of the momentum transfer is generally larger. Therefore, a switch allows 
to choose between BH mass or inverse of the Schwarzschild radius as definition of momentum transfer. 

The part of cm. energy of the neutrino-nucleon system which is not trapped or lost in gravitational radiation at 
formation is attributed to the nucleon remnant. For sake of simplicity, only neutrino-proton collisions are implemented 
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in the MC code. A neutrino-neutron collision does not produce significant statistical differences in the nucleon 
remnant compared to a neutrino-proton collision. The proton remnant is successfully fragmented according to QCD 
in niesons/baryons (see, for example, Ref. |0|) and then is decayed with the PYTHIA generator along with the 
quanta created in the BH evaporation process. 



C. BH evaporation 



Due to the lack of results for the balding and spin-down phases described above, energy losses in these stages are 
assumed to be either negligible or included in the energy loss during formation. This is a reasonable assumption 
since the trapped-surface model likely overestimates the actual energy loss. Balding and spin-down effects are also 
not expected to produce detectable differences in BH air showers, given experimental uncertainties and statistical 
fluctuations. Nevertheless, keeping an open mind, we plan to include balding and spin-down effects in updated 
versions of the code, as soon as theoretical results become available. 

A similar conservative approach is used in the Hawking phase, where only thermally-averaged greybody factors in 
four dimensions are implemented in the MC. This is justified by consistency reasons in the code (the full greybody 
factors for all fields are not known). As the SM fields are emitted on the brane, and given the observational un- 
certainties, the difference between thermally-averaged and exact greybody factors is not expected to be detectable. 
The particle content at trans-Planckian energies is assumed to be the minimal SU{3) x SU{2) x U{1) SM with three 
families and a single Higgs boson. The degrees of freedom Ci and the thermally-averaged greybody factors Ti are 
listed in Table I. The decay multiplicities per species Ni are assigned according to the prescription of Ref. [s^ : 



E,c,r,/,(3) 



(3) 



where fi{m) = 1 or 1 — 2-'^ for bosons or fermions, and the total multiplicity N is 



N 



30C(3)^ E.c»r,/,(3) 



(4) 



where S is the initial entropy of the BH. 

The presence of a minimum length affects the BH evolution in the Hawking phase. If no minimum length is present, 
the MC evaporates the BH accordin g to the Hawking theory. Alternatively, the BH evolution proceeds according to 
the modified thermodynamics of Ref. [i^jyj. In both cases the evaporation ends when the BH reaches the mass Qmin- 
This is set in units of M^, {Mmi) if the minimum length is zero (nonzero). Note that the BH minimum formation mass 
Mmin and the endpoint of Hawking evaporation Qmin are independent parameters. Four-momentum is conserved at 



quarks 

charged leptons 

neutrinos 

photon 

EW bosons 

gluons 

Higgs 

graviton 



72 
12 
6 
2 
9 
16 
1 

(n + 4)(n- 



l)/2 



0.6685 
0.6685 
0.6685 
0.2404 
0.2404 
0.2404 
1 

0.0275 



TABLE I: Degrees of freedom a and thermally-greybody factors Ti for the SM fields. The graviton is assumed to propagate in 
all dimensions. 



each step in the evaporation process by taking into account the recoil of the BH on the brane due to the emission of 
the Hawking quanta. The initial energy of the BH is distributed democratically among all the Hawking quanta with 
a random smearing of ±10%. This smearing factor is introduced on a purely phenomenological basis to take into 
account quantum uncertainties in the emission of each quantum. 
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D. BH final decay 

The MC code allows for two different choices of final BH decay: Final explosion in a number Up of quanta or BH 
remnant. If Up ~ 0, the BH settles down to a remnant with mass Qrniri' If i^p — 1. . . 18, the BH decays in a number 
Up of quanta by a n-body process with total cm. energy equal to Qmin- 

A switch controls conservation of EM charge in the decay process (Hawking evaporation + final decay). The purpose 
of this switch is to allow for the existence of a charged BH remnant. If the EM charge is not conserved and Up = 0, 
the BH remnant carries a charge Qr, where 1 — Qb. = Qh + Qn is the sum of the EM charge of the Hawking quanta 
plus the charge of the nucleon remnant. If the EM charge is conserved and Up — Q, the BH remnant is assumed to be 
electrically neutral, i.e. Qr — 0: The absolute value of the total charge in the Hawking quanta is \Qh\ < 2e/3 and 
Qat = 1 — Qh- This is justified from the fact that the BH charge should have been shed earlier in the evaporation 
process. (See, however, Ref. li^l for a different viewpoint.) It should be stressed that the air shower phenomenology 
of a charge remnant is not known and it is not clear how to track it in the atmosphere in a meaningful way. 



IV. BASICS OF NEUTRINO AIR SHOWERS 

This section presents the essentials of the theory and phenomenology of UHECR neutrino air showers. 



A. Physics of neutrino air sliowers 

UHECRs are believed to be a composite of protons and heavier nuclei. Ultrahigh-energy neutrinos are created 
as these UHECRs interact with the cosmic microwave background through photopion production (protons) or the 
infrared background (iron nuclei) . A cutoff in the energy spectrum is expected at the threshold energy of the photopion 
production, known as the Greisen-Zatsepin-Kuzmin (GZK) cutoff. These neutrinos are called cosmogenic or GZK 
neutrinos. Cosmogenic neutrinos are almost "guaranteed" to exist, though they have not been observed yet. They are 
the most likely source of neutrinos to produce BHs in the atmosphere. The peak of the cosmogenic neutrino fiux is 
around 10^^~^* eV (cm. energy ~ 10 - 50 TeV). The fiux depends on a number of parameters of the UHECR source 
such as spatial distribution, injection spectrum, abundance, maximum energy and cosmological evolution. These 
factors can affect the flux even by a couple of orders of magnitude. 

The depth of the first interaction point, Xq, depends on the total cross section of the process. The column depth of 
Earth's atmosphere in the horizontal direction is 3.6 x 10^ g cm^^. The interaction length of a neutrino with energy 
El, = l(f TeV is \cc — 1-1 x 10^ g cm~^ for charged current (CC) interactions. The largest possible cross sections of 
BH events give shorter interaction lengths, but still larger than the column depth of Earth's atmosphere. Therefore, 
neutrinos interacting in the atmosphere can induce air showers at any Xq. In contrast, SM hadronic interactions have 
large cross sections with Xq high in the atmosphere. Considering deeply penetrating horizontal air showers effectively 
filters out SM hadronic air showers, while giving the most likely chance of interaction. The background for detecting 
BH air showers is limited to SM neutrino air showers. 

The SM interaction channels are the CC and neutral current (NC) for all three flavors. The energy of the leading 
lepton in the flnal state is given by {1 — y)E^, where the mean value of the inelasticity y is w 0.2. The leading lepton 
in the NC interaction is a neutrino that does not contribute to the air shower and the cross section is lower than 
the cross section of the CC interaction, (Jmc{E) ~ QAacc{E). The Vp^-GG produces a high energy /i that does not 
decay before reaching ground. The ^^^.-CC produces a r that also generally does not decay before reaching ground. 
Therefore, the most relevant background for BH air showers is the I'e-CC channel. 



B. Air siiower detection 

Extensive air showers can be detected with fluorescence telescopes and ground arrays. Fluorescence telescopes 
observe the fluorescence light produced by the interaction of atmospheric nitrogen molecules with the EM component 
of the developing air shower. The fluorescence method pioneered by the Fly's Eye detector 41] and currently operated 
by HiRes 42] and the Pierre Auger Observatory (PAO) ^3] is able to reconstruct the longitudinal development of the 
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(mainly) e~^e~ component of the air shower. This technique provides a good estimate of the energy of the primary 
particle that initiates the air shower, since most of the energy of the air shower goes into the observable EM channel. 
This method also enables reconstruction of the shower maximum Xm, i-e. the depth at which the cascade contains 
the maximum number of e+e" pairs. The quantity Xm — Xq is sensitive to the type of primary particle, its energy, 
and the kind of interaction initiating the cascade. The duty cycle is approximately 10%, as clear moonless nights are 
required. 

Ground arrays record the "footprint" of the air shower. Various methods are used to detect charged particles on 
the ground. Some examples are plastic scintillators [i^. l45l| and water Cerenkov tanks [isL l46l| . These detectors are 
spread over areas ranging from a few km'^ to a few thousand km^. For example, the fully functioning southern PAO 
will have 1600 water Cerenkov tanks each with surface area of 10 m^, covering 3000 km^. Arrival time, composition, 
and pattern of the ground signals are used to analyze the properties of the air shower. Ground arrays can be operated 
full time. 

The best method of BH air shower detection is a combination of fluorescence and ground detectors, such as the PAO. 
In view of this, we study both the air shower longitudinal development and the muon content at ground level. For 
inclined air showers, the geomagnetic field affects the distribution of particles on the ground, which is very sensitive 
to the zenith angle. Therefore, we simply count the number of particles rather than study their distribution. 



V. AIR SHOWER GENERATOR 

The BH generator output consists of a list of elementary SM quanta which are decayed with PYTHIA. The 
secondaries of the BH decay (PYTHIA output) are boosted to the laboratory frame and injected in the air shower 
generator AIRES to obtain the air shower. In this section we describe in detail the air shower part of the simulation. 

The AIRES code requires to input primary cosmic ray properties and simulation conditions. The required physical 
parameters are: 

1. Energy of the primary cosmic ray; 

2. Zenith angle of the primary cosmic ray; 

3. Azimuth angle of the primary cosmic ray; 

4. Total number of air shower simulations; 

5. Starting point of the air showers; 

6. Number of observing levels for the longitudinal air shower development; 

7. Observation site to determine geomagnetic field and ground altitude; 

8. Thinning level; 

9. Thinning weight limitation; 

10. Threshold energies for gamma rays, electrons, muons, mesons, and nucleon; 

11. Threshold energies for (approximately) propagating gamma rays and electrons. 

In our simulations, the zenith angle is set to 70° and Xq is set to an altitude of 10 km, corresponding to a slant depth 
of 780 g cm^^. The BH interaction takes place at the injection altitude. The threshold energies for tracking particles 
in the air showers are 100 keV for gamma rays, electrons and positrons, 1 MeV for muons, 1.5 MeV for mesons, and 
150 McV for nucleons. The geomagnetic field is set to the Pierre Auger Observatory (El Nihuil site). The thinning 
level is 10~^ with weight limitation of 0.2. A more detailed explanation of each parameter and other possible options 
can be found in the AIRES manual. 
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A. BH air showers 

The steps to initiate a BH air shower are: 

1. The BH is decayed in the cm. frame. The unstable quanta are hadronizcd or decayed instantaneously by 
PYTHIA, with the exception of top quarks and r leptons. PYTHIA does not handle top quarks. Therefore, 
they are instantaneously decayed as t ^ bW before being injected in PYTHIA. The r leptons are produced 
directly from BH evaporation, or from hadronization or decay of other particles. Depending on their energy, the 
rs may decay before reaching the ground. In that case, they are decayed with PYTHIA but their secondaries 
are injected into AIRES at different atmospheric depths, according to their boost and free path (see below). 

2. All secondaries from PYTHIA are boosted to the laboratory frame. The particles are tightly beamed due to 
their very high boost. 

3. All secondaries are injected into AIRES as primaries of the air shower. 

Note that neutral pions generated in the hadronization process are immediately decayed by PYTHIA in the cm. 
frame. Their average energy in the laboratory frame is smaller than the critical energy, making them more likely to 
decay than interact. 



B. SM air showers 

The CC and NC SM air shower simulations follow Ref. 0: 

1. The differential cross section is integrated over the fraction of the total nucleon momentum carried by the parton 
for all possible values of y. 

2. y is sampled from the previous distribution. The energy of the leading lepton is (1 — y)E^. 

3. The leading lepton of the CC interaction is injected into AIRES. The leading neutrino of the NC interaction is 
not observable and is not injected. 

4. The hadronic part of CC and NC interactions are hadronized with PYTHIA in the cm. frame. The resulting 
particles are boosted back into the laboratory frame and injected into AIRES. 

The T lepton in the i/t--CC interaction is treated separately as in the BH air showers. The r is decayed with PYTHIA 
and its secondaries are injected at the corresponding atmospheric depth. 



VI. SIMULATION RESULTS: SM VS. BH AIR SHOWER DETECTION IN FLUORESCENCE 

DETECTORS 

Simulations show that the characteristics of BH air showers and SM v^-CC air showers with identical first interaction 
point Xq are quite different. However, these differences can vanish with a suitable shift of Xq for either of the two 
processes. As Xq is not a fixed parameter for these interactions, the properties of BH air showers are of limited 
practical use for detection purposes. The characteristics of BH air showers and SM i^e-CC air showers are summarized 
in the following table: 





BH air showers 


SM j^e-CC air showers 


Muon content 


High 


Low 


Development 


Quick 


Slow 


Peak fluctuations 


Small 


Large 


Average total energy 


Varying 


Stable 
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The BH air showers are similar to hadronic air showers. The Ve-GG air showers are comparable to air showers 
generated by photons. The hadronic nature of the BH air showers is due to the prevalence of hadronic channels 
in the Hawking evaporation phase. Their rapidity is due to the large number of hadrons initiating the air showers. 
On the contrary, the main interaction channels in the SM events are pair production and bremsstrahlung. These 
processes produce a smaller number of secondaries than a hadronic interaction. The Landau-Pomcranchuk-Migdal 
effect also contributes in slowing the air shower for primary energies > 10'^ TeV 0. The BH air showers exhibit 
smaller fluctuations in than the Vg-GG air showers. Although the BH mass varies from shower to shower, their 
development is more or less stable because of the large number of BH secondaries. The BH air showers can be viewed 
as a superposition of many air showers with less energy. The larger fluctuations from shower to shower in the SM 
process are due to the fluctuations in the energy (1 — y)E^ carried by the leading lepton. The large variations in the 
total energy of the BH air showers are due to the presence of invisible channels (gravitons, neutrinos, and non-decaying 
Ts). For example, the proton remnant may fragment into a top quark, which decays as t — > 6 + W. If the W decay 
mode is leptonic, a consistent part of the initial proton remnant energy may be carried away by a neutrino. On the 
contrary, most of the cm. energy in the i^e-CC interaction is observable. This leads to a stable air shower total energy. 

The differences between SM and BH air showers can be quantified by choosing a benchmark model for the BH 
process and comparing this model to the SM process. This method also allows to differentiate the effects of various 
parameters and theoretical models from the stable characteristics of the BH air showers. With some guidance from 
the theory, a reasonable choice is: 



Parameter 


Benchmark value 


Planck mass (M*) 


1 TeV 


Number of extra dimensions (n) 


6 


Formation model 


Black disk 


Minimum- allowed mass (Mmin) 


2 


Quantum threshold {Qmin) 


1 


Final hard quanta {rip) 


2 


Momentum transfer (Q) 


i?-i(M) 


EM charge conservation 


YES 


Minimum spacetime length [Im) 






The left panel of Fig. ^ shows the v^-GG air showers and the benchmark BH air showers (50 runs each, neutrino 
primary energy E^, = 10^ TeV). The difference in the shower maxima is AXm = 214 ± 12 g cm^^. Although the 
showers appear to be quite distinct, this difference is a consequence of the same choice of Xq for both BH and SM 
air showers. The right panel of the figure shows i/^-GG air showers shifted so Xm{BH) ~ Xm{ve-GG). Therefore, BH 
and SM air showers can only be distinguished when Xq — Xm is clearly measured. Since present detectors cannot 
measure A"o, v^-GG air showers and BH air showers cannot be discriminated on an event-by-event basis 0. We will 
see in Sect. IVlJTl that this conclusion does not substantially change if the BH parameters are varied. 



VII. SIMULATION RESULTS: SM VS. BH AIR SHOWER DETECTION IN GROUND AND HYBRID 

DETECTORS 



Discrimination of BH air showers and v^-GG air showers can be improved by the use of ground arrays. The best 
possible scenario for BH air shower detection is a technique that combines air fluorescence telescopes and a ground 
array. Since fluorescence telescopes are able to measure Xm accurately, a good air shower discriminator is to fix Xq 
and count the number of particles at various distances from X.„i. This is equivalent to fixing the detection level and 
varying the air shower first interaction point. The fluctuation due to the change of Xq is negligible compared to 
fluctuations arising from other uncertainties. Figure |21 shows the number of muons at different atmospheric depths vs. 
the number of electrons at the air shower maximum for BH air showers (benchmark model) and i/g-CC air showers. 
The muons are measured from the ground array and the electrons are measured from the fluorescence telescopes. The 
BH air showers are characterized by a higher muon content than the Vf.-GG air showers. Although the separation is 
not large enough to distinguish the air showers on an event-by-event basis, discrimination of BH and SM events is 
possible with enough statistics. The number of muons also depends on the detection level. Since different AXs can 
be seen as different initial interaction points, the air showers in the lower right (upper left) panel of Fig. [3 can be 
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FIG. 1: Number of e e~ vs. slant depth for the longitudinal development of 50 air showers with = 10 TeV. The BH air 
showers for the benchmark model (black solid lines) and the Ve-CC air showers (red dashed lines) are shown. The air shower 
maxima are Xm = 1566 ± 6 g cm~^ for the BH benchmark model and Xm ~ 1780 ± 9 g cm~^ for the v^-CC, respectively. 
The difference in the air shower maxima is AX™ = 214 ± 12 g cm~^. The left panel has both air showers with identical first 
interaction point, Xo(ve-CG) = Xo{BH) — 780 g cm~^. The right panel shows the same air showers with a shift in Xoive-CC) 
such that Xmii^e-CC) ~ X^{BH). 
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FIG. 2: Number of ^i^ fJ." at various atmospheric depths Xm + AX vs. the number of e^e~ at Xm for 50 benchmark model BH 
air showers (black filled circles) and 50 Ve-CC air showers (red empty circles). The energy of the primary neutrino is E^, — 10^ 
TeV. The observation depth of the muons increases from left to right panel and from top to bottom panel. 
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FIG. 3: Number of j-i'^ fJ-~ vs. number of e^e~&t various depths Xm + AX. The comparison is between 50 BH benchmark 
model air showers (black filled circles) and 50 Ve-CC air showers (red empty circles) with E^, — 10^ TeV.The observation depth 
of the muons and electrons increases from left to right panel and from top to bottom panel. 



understood as starting higher (lower) in the atmosphere than the air showers in the other panels. The number of 
muons in the air showers decreases if the first interaction point is higher in the atmosphere. 

If only ground detection is possible, a good BH vs. SM discriminator is the number of electrons and muons at 
various atmospheric depths X,„ + AX (Fig. As in the hybrid detection scenario, although the BH air showers 
show higher muon content than the Vf.-CC air showers, their discrimination requires large statistics. The number of 
muons depend on the observation level and decreases as AX increases. In absence of an air fluorescence telescope to 
accurately measure X^, alternative techniques must be used to reconstruct the air shower maximum. 



VIII. SIMULATION RESULTS: EFFECT OF BH PARAMETERS 



The effects of various BH models can be studied by varying the parameters described in Sect. Illll and comparing 
the simulations to the benchmark model. The air shower longitudinal developments for different choices are shown in 
Figs. I4I8I and summarized in Table II. 

The depth of the air shower maximum is not significantly affected by changes in the BH parameters or the theoretical 
model. Fluctuations in the air showers are large; all X„i values lie within one standard deviation from each other, with 
the possible exception of the simulations with nonzero minimum length. We conclude that the main characteristics 
of BH air showers described in the previous section are robust. This result can be qualitatively explained by noticing 
that most of the cm. energy of the neutrino-nucleon collision is not trapped in the BH. Therefore, different choices of 
BH parameters do not produce large observable effects in the air shower development. The main factors determining 
the BH evolution (BH mass distribution, energy and spectrum of emitted quanta) are difficult to disentangle because 
their variations do not act to coherently increase or decrease the shower maximum. For instance, increasing the 
minimum-allowed BH mass from M„vin — 2 TeV to Mmin = 10 TeV increases the average BH mass in the air 
showers. This leads to a larger number of quanta. However, this property does not translate into a faster air shower 



13 



Difference from benchmark 


Xm ± rms 


error 


None (benchmark) 


1566 


± 


6 


M* = 3 TeV 


1545 


± 


6 


n = 3 


1551 


± 


5 


Mrr^iu = 10 TeV 


1546 


± 


6 


Qmin — 2 TeV 


1559 


± 


7 


Neutral remnant 


1549 


± 


6 


Charged remnant 


1564 


± 


6 


YN model 


1547 


± 


6 


YR model 


1547 


± 


5 


Imin ~ 2.5 


1519 


± 


4 



TABLE II: Shower maximum Xm and rms error for 50 BH air showers with different physical parameters and models. The 
primary neutrino energy is E^, — 10^ TeV. The first row gives Xm for the benchmark model of Sect. IVll 
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FIG. 4: Longitudinal development of 50 air showers for the BH benchmark model (black solid curves) vs. two different choices 
of BH parameters (red dashed curves). The energy of the primary neutrino is E^, = 10^ TeV. The left panel shows the 
difference between the benchmark model (M* = 1 TeV) and = 3 TeV. The average difference in the air shower maxima 
is AXm = 21 ± 9 g cm~^. The right panel shows the benchmark case (n = 6) and n = 3. The average difference in Xm is 
AX„ = 15±10 gcm-^ 

development; the average energy per quanta is smaller, and the two effects compensate each other. Two interesting 
facts are worth observing. Firstly, the benchmark case has the largest cross section and the largest X„i . This result is 
mainly due to the choice of a relatively small fundamental Planck constant. Adding quantum effects or graviton loss 
at formation seems to decrease X^n slightly. Secondly, the presence of a minimum length may possibly be the only 
BH physical signature distinguishable from the black disk model. However, in our simulations the choice of Imin has 
been purposedly fine-tuned to the maximum-allowed value that allows BH formation with primary neutrino energy 
Eu = 10^ TeV. Relaxing this choice leads to values of Xm closer to the black disk result. 



IX. SIMULATION RESULTS: BACKGROUND-FREE SIGNAL FROM r DECAY 

It has been suggested that r leptons produced during BH events, cither directly or through decay of other particles, 
can produce an observable signal similar to the double bang produced in Vt-CC interactions 10.]. The mechanism is 
the following. The BH evaporation initiates a first air shower, while a second air shower is initiated by the t decay at a 
lower altitude. If the second bang is large enough to be observed, the double peak feature provides a background-free 
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FIG. 5: Longitudinal development of 50 air showers for the BH benchmark model (black solid curves) vs. two different choices 
of BH parameters (red dashed curves). The energy of the primary neutrino is = 10^ TeV. The left panel shows the 
difference between the benchmark model {Mmin = 2Mi,) and Mmin = lOM*. The average difference in the air shower maxima 
is AXm = 20 ± 10 g cm~^. The right panel shows the benchmark case (Qmin = M*) and Qmin = 2Mi. The average difference 
in X„ is AX™, = 7 ± 10 g cm"^ 
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FIG. 6: Longitudinal development of 50 air showers for the BH benchmark model (black solid curves) vs. two different choices 
of BH parameters (red dashed curves) . The energy of the primary neutrino is Ei, = 10^ TeV. The left panel shows the difference 
between the benchmark model (final decay in 2 quanta) and BH evolution with final electrically neutral remnant. The average 
difference in the air shower maxima is AXm = 17 ± 10 g cm"'^. The right panel shows the benchmark case and BH evolution 
with final electrically charged remnant. The average difference in Xm is AXm = 2 ± 10 g cm~^. 



signature independent of the first interaction point. 

The T lepton in the Vr-CC interaction carries on average about 80% of the total cm. energy. In a BH event, the 
fraction of the total BH mass going into the r second bang can be estimated from Eq. © to be roughly 

lr,^, -^%- (5) 



Therefore, a r from BH evaporation carries on average a few percent of the total BH mass. Since the energy trapped 
in the BH accounts generally only for a small fraction of the initial cm. energy, the total energy in the t channel is 
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FIG. 7: Longitudinal development of 50 air showers for the BH benchmark model (black solid curves) vs. two different choices 
of BH[ parameters (red dashed curves). The energy of the primary neutrino is Ei, — 10^ TeV. The left panel shows the difference 
between the benchmark model (black disk) and the YN graviton loss model. The average difference in the air shower maxima 
is AXm = 19 ± 11 g cm~^. The right panel shows the benchmark case and the improved YR graviton loss model. The average 
difference in Xm is AXm = 19 ± 8 g cm~^. There is virtually no difference between the YN and YR models. 
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FIG. 8: Longitudinal development of 50 air showers for the BH benchmark model (black solid curves) and a BH evolution 
model with minimum length Imin = 2q: M^^ (red dashed curves). The energy of the primary neutrino is Ei, — 10^ TeV. The 
average difference in the air shower maxima is AXm = 47 ± 13 g cm~^. 



on average less than 1% of the total air shower energy. This gives the t a higher chance to decay before reaching 
ground, but a smaller energy deposit in the second bang, making the latter harder to detect. Moreover, rs are not 
produced every time a BH is formed and the probability to observe an air shower with a r decaying before reaching 
ground is relatively small. 

The double bang signature can be studied by selecting only air showers containing at least one r decaying in air. 
We simulated these air showers using the benchmark parameters, zenith angle 70° and Xq = 160 g cm~^ (altitude 
of the first interaction point = 20 km). The higher altitude gives a larger separation between the two bangs. Two 
possible scenarios were considered: i) rs decaying at any depth in the atmosphere, and ii) at least one r decaying at 
an altitude X > 0.75{Xg — Xq), where Xg is the slant depth of the ground ("low altitude rs"). Case ii) represents 
the best possible scenario for double bang detection, as the second bang is expected to occur close to the detector. 
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FIG. 9: Longitudinal development of 50 BH air showers containing at least one r decaying at an altitude X > 0.7^{Xg — Xo). 
The air showers are started at an altitude of 20 km, corresponding to a slant depth of Xq — 160 g cm~'^. The left panel shows 
the profile of all the particles in the air shower. The right panel shows only the profile of the secondary particles from the r 
decay. Note that the y sxis of the right panel is magnified 10 times w.r.t the left panel. The peaks between 500 g cm~^ and 
1500 g cm~^ are originated by multiple r events with one r decaying high in the atmosphere. 

The longitudinal profiles of the low altitude r air showers are shown in the left panel of Fig. |^ The second bang is 
visible around 2250 g cm~^. The longitudinal development of the r channel contribution to the air shower is plotted 
in the right panel. The peaks between 500 g cm~^ and 1500 g cm~^ are from events with multiple rs, where at least 
one of these ts decays at low altitude. As is expected, the r component contributes only a minimal fraction to the 
overall air shower energy. Although the second bang is in principle detectable, it is within the fluctuation of the first 
bang's tail; current and next generation detectors cannot discern a few percent feature. Moreover, existing detectors 
such as the PAO have a limited field of view and cannot track the full profile. 



X. CONCLUSIONS 

We simulated extensive air showers initiated by TeV-scale BH events produced from neutrino interaction in Earth's 
atmosphere. These simulations were performed with an advanced fortran MC code 16] that includes most of the 
theoretical results of the recent literature. The BH air showers were compared to SM air showers and different models 
of BH formation and evolution were investigated. We also studied exotic signatures of BH events, such as the "r 
double bang" signature. Our goal was to test various proposals of BH vs. SM air shower discrimination and look for 
new ways of differentiating models of BH formation and evolution. 

Our results show that the main features of BH air showers are largely independent of the details of BH formation 
and evolution. Statistical fluctuations and limitations in detection techniques hinder the discrimination of alternative 
theoretical or phenomenological models. No difference between the black disk model and alternative models of BH 
formation, or between BHs with final explosive decay and stable remnant, can be detected with current UHECR 
observatories. Distinguishing SM and BH air showers with hybrid detectors is possible if enough statistics is gathered. 
The most promising way is to measure the air shower maximum with a fluorescence telescope and count muons at 
ground. The double bang signature cannot be observed by any realistic detector at the present stage. These results 
imply that the theory of TeV-scale BHs in UHECRs is robust, but BH air showers are hardly to probe details of "new 
physics" in the near future. 
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